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We formulate a theory of collective phase description for oscillatory convection in Hele-Shaw 
cells. The theory can be considered a phase reduction method for limit-cycle solutions in infinite- 
dimensional dynamical systems, namely, time-periodic solutions to partial differential equations 
representing oscillatory convection. We derive the phase sensitivity function, which quantifies the 
, phase response of an oscillatory convection to weak perturbations applied at each spatial point. 

■ Using the phase sensitivity function, we analyze phase synchronization between two weakly coupled 

oscillatory convections. 
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INTRODUCTION 



Self-sustained oscillations and synchronization phenomena are ubiquitous in nonlinear systems, which are typically 
described by ordinary differential equations [iHll • Nonlinear systems also possess abundant examples of spatiotcmpo- 
ral pattern formation, which are generally described by partial differential equations, e.g., reaction-diffusion equations 
and fluid equations 0,01 • Naturally, synchronization phenomena between spatiotemporal patterns (e.g., spatiotem- 
poral chaos) have attracted considerable attention among researchers 

Examples of earlier studies in reaction-diffusion systems include the following. Experimental and theoretical studies 
have been carried out on synchronization between two locally coupled domains of excitable media exhibiting spiral 
■ wave behavior using the photosensitive Belousov-Zhabotinsky reaction [lol|. Synchronized pulses in laterally coupled 
excitable fibers have been numerically analyzed using the spatially one-dimensional FitzHugh-Nagumo equations . 
The dynamics of propagating waves under spatiotemporal forcing have been investi gat ed analytically and numerically 
for a model system undergoing phase separation as well as chemical reactions |12Hl4j . 

There are also examples in fluid systems, which are mainly aimed at understanding synchronization phenomena in 
geophysical fluid dynamics. Numerical studies have been performed on synchronization between spatiotemporal chaos 
in a pair of quasi-two-dimensional channel models [15l] and in a pair of Hele-Shaw cells (Tcj . In addition, experimental 
studies have been carried out on synchronization in both periodic and chaotic regimes in a pair of thermally coupled 
rotating baroclinic annulus systems [13, [3 ■ 

In Refs. [l9l - l22j . we recently formulated a theory of collective phase description for globally coupled noisy dynamical 
elements exhibiting macroscopic oscillations, which reduces the nonlinear Fokker-Planck equation (a partial integro- 
differential equation) to the collective phase equation (an ordinary differential equation). In particular, in Ref. [23], 
we considered the nonlinear Fokker-Planck equation describing globally coupled noisy active rotators, which does 
not possess spatial translational symmetry. The theory can be considered a phase reduction method for limit- 
cycle solutions in infinite-dimensional dynamical swtems. A similar formulation for stable time-periodic solutions to 
reaction-diffusion systems can also be developed [23|. 

In the present paper, as in Ref. p^ . we consider synchronization between Hele-Shaw cells [41*1. Focusing our atten- 
tion on stable time-periodic solutions, we formulate a theory of collective phase description for oscillatory convection 
in Hele-Shaw cells, which enables us to describe the dynamics of an oscillatory convection by a single degree of free- 
dom. The theory can be considered an extension of our phase reduction method for the nonlinear Fokker-Planck 
equation [23| to an equation for oscillatory Hele-Shaw convection. On the basis of our theory, we analyze phase 
synchronization between two weakly coupled oscillatory Hele-Shaw convections. 
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The notion of collective phase originates from the phase of the collective oscillation emerging from coupled indi- 
vidual phase oscillators (T^-Ull- In the present paper, as in Ref. j22|, the collective phase is assigned to temporal 
translational symmetry breaking in partial differential equations. In general, the phase arises not only from tempo- 
ral translational symmetry breaking but also from spatial translational symmetry breaking [3]. In fact, the phase 
dynamics of spatially periodic structures have been extensively developed, which are based on spatial translational 
symmetry breaking [2J-l27| . There is a history of the phase dynamics approach to spatially periodic patterns in 
fluid systems [28l-34[ (see also Refs. [1, [1]). In addition, the so-called interface dynamics or pulse dynamics are also 
essentially based on spatial translational symmetry breaking 35, 36]. In contrast to these studies, our formulation is 
based on temporal translational symmetry breaking. Therefore, the formulation is applicable to oscillatory Hele-Shaw 
convection, although this system does not possess spatial translational symmetry. 

The organization of the present paper is as follows. In Sec. [TTl we formulate a theory of collective phase description 
for oscillatory Hele-Shaw convection. In Sec lIIIi we illustrate our theory using numerical simulations of the oscillatory 
convection. Concluding remarks are given in Sec. IIVI 



II. PHASE DESCRIPTION OF OSCILLATORY HELE-SHAW CONVECTION 

In this section, we formulate a theory of collective phase description for oscillatory Hele-Shaw convection. The 
theory can be considered an extension of our phase reduction method for the nonlinear Fokker-Planck equation [22j | 
to an equation for oscillatory convection. 

A. Dimensionless form of governing equations 



We consider the dynamics of the temperature field in the Hele-Shaw cell to be described by the following dimen- 
sionless form (see Ref. [16] and also references therein): 



d_ 

dt 



Tix,y,t)^V^T + Ji,p,T), 



where the Laplacian and Jacobian are respectively given by 



dT di: dT 



v^r = ( — — \- ^— I r. 



dx dy dy dx 

The stream function il){x^y,t) is determined from the temperature field T(x,y,t) as follows 

dT 



V'^{x,y,t) 



-Ra 



dx ' 



where the Rayleigh number is denoted by Ra. The stream function also provides the fluid velocity fleld, i.e. 

, , f d-d) d'ilj\ 



(1) 

(2) 
(3) 

(4) 
(5) 



The system is defined in the unit square, i.e., x G [0, 1] and y G [0, 1]. The boundary conditions for the temperature 
field T{x, y, t) are given by 



dT{x,y,t) 



dx 

T{x,y,t) 



dT{x,y,t) 



x=0 
= 1, 



dx 
T{x,y,t) 



-0, 



= 0, 



(6) 
(7) 



where the temperature at the bottom {y = 0) is higher than at the top {y — 1). The stream function ^l^{x, y, t) satisfies 
the Dirichlet zero boundary condition on both x and y, i.e.. 



ip{x,y,t) 

ip{x,y,t) 



x=0 



'il^{x,y,t) 
'ilj{x,y,t) 



x=l 



0, 



0. 



(8) 
(9) 
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B. Variational components of temperature fields 

To simplify the boundary condition ([7]), we consider the variational component of the temperature field as follows 



T{x,y,t) {l-y)+X{x,y,t). 
Inserting Eq. (fTO| into Eq. ^ and Eq. we derive the following equations: 

and 

dX 

V^i^{x,y,t) = -Ra— . 
Applying Eq. pO)) to Eqs. ([6])([71), we obtain the following boundary conditions: 



dX{x,y,t) 



dx 



dX{x,y,t) 



x=0 



dx 



= 0, 



x=l 



X{x,y,t) =X{x,y,t) 
v=o 



= 0. 



(10) 

(11) 
(12) 

(13) 
(14) 



Namely, the field X{x,y,t) satisfies the Neumann zero boundary condition on x and the Dirichlet zero boundary 
condition on y. 

In the manipulation below, it should also be noted that Eq. (jl2p can be rewritten in the following form: 



'^{x,y,t)^J dx' J dy' G{x,y,x',y')-^X{x',y',t), 

where the Green function G(x,y,x\y') is the solution to 

V^G(x, ?/, a;', y') = ~Ka5{x - a;') 5{y - y'), 



(15) 



(16) 



under the Dirichlet zero boundary condition on both x and y. In the following two subsections, we analyze the 
dynamical equation (fTT|) with Eq. (fT2|) or Eq. (|15p . under the boundary conditions given by Eqs. ([T3l)(|14|) and 
Eqs. dSD®. 



C. Time-periodic solution and its Floquet-type system 

The stable time-periodic solution to Eq. (|lip . which represents oscillatory convection in the Hele-Shaw cell, can be 
described by 



X{x,y,t)=XQ{x,yMt)). 0(1) ^ n, 



(17) 



where Q and fl are the collective phase and collective frequency, respectively [42]. Inserting Eq. ^T7\ into Eq. (fTTj) 
and Eq. ([T^ . we find that Xo{x, y, &) satisfies the following equation: 



n^Xoix, y, e) = v^Xo + J(v^o, ^o) - 

o<3 ax 



where 



VVo(2;,y,e) = -Ra 



dXp 
dx 



(18) 



(19) 



Let u{x, y, O, t) represent a small disturbance to the time-periodic solution and consider a slightly perturbed solution 

Xix,y,t)=Xo{x,y,eit))+u{x,y,e{t),t). (20) 
Equation is then linearized in u{x,y,Q,t), i.e., 

d 



dt 



u{x, y, 6, t) = C{x, y, Q)u{x, y, 6, t). 



(21) 
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Here, the linear operator C{x, y, Q) is explicitly given by 



C{x,y,Q)u{x,y,&) = 



L{x,y,e)-n 



de 



u{x,y,Q), 



where 



L(a;,y, Q)u{x,y, 6) = V^w + J{il^o,u) + J(V'„,Xo) 



dx 



(22) 



(23) 



Similarly to the time-periodic solution X(){x^y,Q), the function u{x,y,0) satisfies the Neumann zero boundary 
condition on x and the Dirichlet zero boundary condition on y. In Eq. (I23p . the function "ipuix, y, 6) is the solution to 



V'^ipu{x,y,Q) = -Ra.— , 



(24) 



under the Dirichlet zero boundary condition on both x and y. Note that C{x^ y, Q) is time-periodic through 0, and 
therefore, Eq. (j2ip is a Floquet-type system with a periodic linear operator. 
Defining the inner product of two functions as 



u*ix,y,e), u{x,y,e) 



1 

2^ 



dO dx dy u*{x,y,Q)u{x,y,&), 



we introduce an adjoint operator of the linear operator C{x, y, O) by 

u*{x,y,Q), C{x,y,Q)u{x,y,Q) = C*{x,y,<d)u*{x,y,<d), u{x,y,Q) 
By partial integration, the adjoint operator C*{x,y, O) is explicitly given by 



where 



C*ix,y,e)u*{x,y,e) = 



L*(x,y,e)u*(x,y,e) = VV + ^ 



L*{x,y,e) + n 



d_ 

de 



u*ix,y,e), 



. dtpo 
dy 



d_ 

dy 



dx 



d_ 

dx 



(25) 



(26) 



(27) 



(28) 



Similarly to u{x^ y, 0), the function u*{x, y, Q) satisfies the Neumann zero boundary condition on x and the Dirichlet 
zero boundary condition on y. In Eq. (|28p . the two functions, i.e., ip* ^{x, y, O) and ip* y{x, y, 0), are the solutions to 



V2V';,(a;,y,e) = -Ra^ 
VV;,(a;,y,e)--Ra^ 



dx 



(29) 
(30) 



under the Dirichlet zero boundary condition on both x and y, respectively. Details of the derivation of the adjoint 
operator C*{x, y, Q) are given in App. [Xj 



D. Floquet zero eigenfunctions 



In the calculation below, we use the Floquet eigenfunctions associated with the zero eigenvalue, i.e., 

d 



C{x,y,e)Uo{x,y,e) 
£*{x,y,Q)U*{x,y,e)^ 



Lix,y,e)~n 

L*{x,y,e) + n 



de 

d_ 

de 



f/o(x,2/,e) = 0, 
U*ix,y,e)^0. 



Note that the right zero eigenfunction U(){x, y, e) can be chosen as 

d 

Uo{x,y,e) = —Xo{x,y,e), 



(31) 
(32) 

(33) 



5 



which is confirmed by differentiating Eq. ([T8)) with respect to Q. Using the inner product (1251) with the right zero 
eigenfunction (1551) . the left zero eigenfunction J7q(x, y, 0) is normahzed as 



US{x,y,&),Uo{x,y,e)]=^ [ dQ ( dx f dyUS{x,y,&)Uo{x,y,e)^l. 
-11 Jo Jo Jo 

Here, we should note that the following equation holds (see also Refs. [22l. |37|): 



_d_ 



dx / dyUo{x,y,&)Uoix,y,Q) 
Jo 

1 /•! r 5 Q 

USi^, y. ^)^Uo{x, y, 8) + Uo{x, y, e) — US{x, y, 9) 



dx I dy 
Jo 
1 fi 



= ^ / dx I dy 
= 0. 



C/o* (x, y, e)L(x, y, e)[/o(x, y, 9) - [/o(x, y, 9)L*(x, y, 9)l/o* (x, y, 9) 



Therefore, the following normalization condition is satisfied independently for each 9, i.e., 



dx I dyU^{x,y,Q)Uo{x,y,Q) ^ I. 



(34) 



(35) 



(36) 



In the following two subsections, using the time-periodic solution and its Floquet zero eigenf unctions, we formulate a 
theory of collective phase description for oscillatory Hele-Shaw convection. 



E. Oscillatory convection with weak perturbations 



In this subsection, we consider oscillatory Hele-Shaw convection with weak perturbation to the temperature field 
that is described by the following equation: 



^X(x, y, t) = y^X + Jiij, X)~^+ eCix, y, t). 



(37) 



The weak perturbation is denoted by e£,{x,y,t). Using the idea of phase reduction 0, we can derive a phase equa- 
tion from the perturbed equation p7p . Namely, we project the dynamics of the perturbed equation (j37p onto the 
unperturbed solution as 



Q{t) = 1^ dx dy U*{x, y, e)^X(x, y, t) 
c^^n + e dx dyUQ{x,y,Q)^{x,y,t), 



(38) 



Jo Jo 

where we approximated X{x,y,t) by the unperturbed solution Xo(x, y, 9) and used the fact that 

dx f dyf/o*(x,y,9)^Xo(x,y,9) = r! / dx f dy [/*(x, y, 9)[/o(x, y, 9) = f!. (39) 
Jo '^t Jq 

Therefore, the phase equation describing the oscillatory Hele-Shaw convection with weak perturbation is obtained in 
the following form: 



e{t) = n + ef dx [ dyZ{x,y,e)^{x,y,t), 
Jo Jo 

where the phase sensitivity function is defined as 

Z{x,y,e)^U*{x,y,e). 



(40) 



(41) 



It should be noted that Eq. (PO]) corresponds to a phase equation that is derived for a perturbed limit-cycle oscillator 
described by a finite-dimensional dynamical system (see Refs. [l|, 0, 113, H^). However, refiecting the aspect of an 
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infinite-dimensional dynamical system, the phase sensitivity function Z{x^ y, Q) of the oscillatory Hele-Shaw convection 
possesses infinitely many components that are continuously parameterized by the two variables, i.e., x and y. 

Here, we briefly describe a numerical method for obtaining the phase sensitivity function. From Eq. (|32|) , the phase 
sensitivity function Z(x, y, Q) satisfies 

^-^Zix, y, 9) = -L*ix, y, e)Z{x, y, 8), (42) 

which can be transformed into 

^Z(x, y, -fls) = L* (.T, y, -ns)Z{x, y, -ils), (43) 

where Q = —fls. A relaxation method using Eq. (I43p is convenient for obtaining the phase sensitivity function (see 
also Refs. [H, [H, [13, SI ) . 



F. Weakly coupled oscillatory convections 



In this subsection, we consider weakly coupled oscillatory Hele-Shaw convections that are described by the following 
equation (Tsj : 

^^X,{x,y,t) = V^X, + J(^.,X,) - ^^+e{Xr - X,), (44) 

for (cr, r) = (1,2) or (2,1), where the stream function for each system is determined by 

dX 

VV.(a^,2/,i) = -Ra^. (45) 

Two identical oscillatory Hele-Shaw convections are mutually coupled through corresponding temperatures at each 
spatial point 4^], where the coupling parameter is denoted by e. Here, we assume that unperturbed oscillatory Hele- 
Shaw convections are stable time-periodic solutions and the coupling between the oscillatory convections is sufficiently 
weak. Under this assumption, as in the preceding subsection, we can obtain a phase equation from Eq. (|44p as follows 



e,{t)=n + ej^ dxj^ dy Zix,y,e,)(^Xoix,y,er) - Xo{x,y,e,)y (46) 
Applying the averaging method ^ to Eq. dUl), we can derive the following phase equation: 

e,{t) = n + fT{e„-er), (47) 

where the phase coupling function is given by 



r(e) = ^^ dxj^ dxj^ dy Z{x,y,X + e)(^Xo{x,y,\) ~ Xn{x,y,X + 0))■ 



m 



The type of the coupling (e.g., in-phase or anti-phase) is determined by the ant i- symmetric component of the phase 
coupling function Q, i.e., 

ra(e) = r(e)-r(-e). (49) 

Finally, we note that the form of Eq. (j47p is the same as that of the phase equation which is derived from weakly 
coupled limit-cycle oscillators described by finite-dimensional dynamical systems (see Ref. [2]). That is, an oscillatory 
convection can be reduced to a phase oscillator, similarly to an ordinary limit-cycle oscillator. 



III. NUMERICAL ANALYSIS OF OSCILLATORY HELE-SHAW CONVECTION 



In this section, using numerical simulations of oscillatory Hele-Shaw convection, we illustrate the theory developed 
in the preceding section. 
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A. Spectral transform and order parameters 

Considering the boundary conditions of X{x, y, t), i.e., Eqs. (fT3)) ([T4)) . we introduce the following spectral transform: 

Hjk{t) = dx dy X{x,y,t) cos{T:jx) sm^nky), (50) 



for j — 0, 1, 2, • • • and k — 1,2, ■ ■ ■ . In visualizing the limit-cycle orbit in the infinite-dimensional state space, we 
project the time-periodic solution XQ{x,y,Q) onto the H11-H22 plane as 

77ii(6) = dx dy Xo{x,y,Q)cos{TTx)sm{TTy), (51) 

H22{Q)^ f dx f dyXo{x,y,e)cos{2Trx)sm{2TTy), (52) 
Jo Jo 

which can be considered a pair of order parameters. 

B. Time-periodic solution and phase sensitivity function 

First, we summarize our numerical simulations in the present paper. We applied the pseudospectral method that 
is composed of the sine expansion with 128 modes for the Dirichlet zero boundary condition and the cosine expansion 
with 128 modes for the Neumann zero boundary condition. The initial values were prepared so that the system 
exhibits a single cellular oscillatory convection. The Rayleigh number was fixed at Ra = 480, and then the collective 
frequency is ~ 622. 

Figure [1] shows the limit-cycle orbit projected on the H11-H22 plane, which was obtained from our numerical 
simulations of the dynamical equation (1111) . Snapshots of the time-periodic solution Xo{x, y, Q) and other associated 
functions, i.e., To{x, y, 8) = (1 — y) -f Xq{x, y, 8), iioix, y, 8), Uq{x, y, 8), and Z{x, y, 8), are shown in Fig. [2] The 
phase sensitivity function Z{x, y, 8) was obtained using the relaxation method for Eq. (|43| . In addition. Fig. [3] shows 
typical shapes of both the time-periodic solution Xo{x, y, 8) and the phase sensitivity function Z{x, y, 8) with respect 
to 8. 

Here, we describe the property of the phase sensitivity function Z{x, y, 8) under this condition. For each 8, similarly 
to the time-periodic solution Xo(x,y,8), the phase sensitivity function Z{x,y,@) is anti-symmetric with respect to 
the center of the system, i.e., 

Xo{-xs, -ys, 8) = -Xo{xs,ys, 8), (53) 
Z(-a;5, -j/5,8) = -Z(a;5,2/5,8), (54) 

where xs ^ x — 1/2 and ys = y — 1/2. Therefore, the phase sensitivity function at the central point is equal to zero, 
i.e., Z{x = 1/2, y = 1/2, 8) = 0. In addition, the spatial integral of the phase sensitivity function also becomes zero, 
i.e., 

dx [ dyZ{x,y,Q) = 0. (55) 



Jo 



Namely, when the weak perturbation is spatially uniform, i.e., ^(x, y, t) = £,{t), the collective phase is neither advanced 
nor delayed by the perturbation. It should also be noted that the phase sensitivity function Z{x,y,Q) is spatially 
localized. That is, the amplitudes of the phase sensitivity function Z{x,y, 8) with respect to 8 in the top-right and 
bottom-left corner regions of the system are much larger than in the other regions. 



C. Phase synchronization between two weakly coupled oscillatory convections 

In this subsection, as in Ref. [l3|, we consider two weakly coupled oscillatory Hele-Shaw convections that are de- 
scribed by the partial differential equation ()44|) . Here, we assume that unperturbed oscillatory Hele-Shaw convections 
are stable time-periodic solutions and the coupling between the oscillatory convections is sufficiently weak. Under 
this assumption, we can theoretically analyze the phase synchronization between oscillatory convections. 
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The anti-symmetric component of the phase couphng function is shown in Fig. IH^a), which was calculated by 
Eqs. As can be seen in Fig. IDJa), the type of the phase coupling function is in-phase (attractive), i.e., 

T'^{0) < and r^(±7r) > 0. Namely, the two oscillatory convections will become in-phase synchronized. 

Figure lU^b) shows the time evolution of the collective phase difference \Qi — 02 1 between the two oscillatory 
convections, which started from an almost anti-phase state with the coupling parameter e = 0.05. The two oscillatory 
convections eventually became in-phase synchronized, namely, Qi = 02- Our direct numerical simulation of Eq. (|44p 
is compared to the theory, i.e., O — era(O). The simulation result agrees perfectly with the theory. 



IV. CONCLUDING REMARKS 



We developed a theory of collective phase description for oscillatory convection in Hele-Shaw cells, by which an 
oscillatory convection can be reduced to a phase oscillator. On the basis of our theory, we analyzed phase synchro- 
nization between two weakly coupled oscillatory Hele-Shaw convections. The key quantity of our theory is the phase 
sensitivity function of the oscillatory Hele-Shaw convection, which quantifies its phase response to weak perturbations 
applied at each spatial point. 

Finally, we remark the broad applicability of our approach, which is not restricted to the oscillatory Hele-Shaw 
convection. If we focus on stable time-periodic solutions, similarly to ordinary differential equations, partial differential 
equations can generally be reduced to phase equations by our approach. There are abundant examples of rhythmic 
phenomena in nature that can be described by partial differential equations, e.g., in geophysical fiuid dynamics. The 
phase description approach has the capability to play a central role in such areas. 
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Appendix A: Derivation of the adjoint operator 



Here, we describe the details of the derivation of the adjoint operator C*{x,y, 6) given in Eqs. (|27|)(f28|) (see also, 
e.g., Refs. [H,!!^ on mathematical terms). From Eqs. ([^(P^ . the linear operator C{x,y,Q) is given by 



C{x,y, Q)u{x,y, 6) 



d^u d'^u d'ipo du dijjo du d^u f dX^ 



dx'^ dy^ dy dx dx dy dx \ dy 



dipu dXp _ ^dv_ 
dy dx dQ 



By partial integration, each term of the inner product |u*(x, y, O), C{x, y, 0)it(a;, y, O)] can be transformed into 



u* — 

' dx'^ 

' dy'^ 
dipo du 



dy dx 
dipQ du 



dx dy\ 

, d^pu f dXp _ _^ 
' dx \ dy 

djpu dXo 



dy dx 
du 



-n- 



de 



1 

1 

2^ 



dQ / dy 



dQ / dx 



du 
dx 

du 
dy 



x=0 



1 



2tt 



1 

2^ 



dy 



1 

2^ 









L 


/ dx 


u 




'o 




p27T 






L * 


f dy 


u* 




Jo 





dy 

dtpo 



x = l 



du* 
dx 

du* 
dy 



y=o. 



d^u* "n 

d^u* ' 
dy'^ ' 



+ 



1 



2tt 



dx 



dx 

f dXo 
V dy 

. dXp 
dx 



x=0 

y=i 



y=o 

y=i 



d_ 

dx 

d_ 

dy 



. dijjQ 
dy 

dy 



2^ 



dx / dy 



y=o 
e=2TT rr 



x=o iL 

d_ 

dy 



_d_ 

dx 

dXo 
dx 



dXo 
dy 



0=0 



du* 



(Al) 

(A2) 
(A3) 
(A4) 
(A5) 
(A6) 
(A7) 
(A8) 
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Using the Green function G{x, y, x' , y') given in Eq. (|16p . the function ipuix, y, 9) given in Eq. (p4|) can be rewritten 
in the following form: 



ipu{x,y,Q)= / dx' / dy' G{x,y,x' ,y') 



^0 



, , du{x',y',e) 



dx' 



In Eqs. (|A6|)(|A7[) . we perform the following manipulations: 

'dXo 



d_ 

dx 



dy 



1 



1 



dQ / dx dy 



^0 

1 



d_ 

dx 



/dXo 
V 9y 



- 1 



-hi '""l ''I '^io ^'^S^^ 







de dx dy dx' / dy' G{x',y',x,y) 



1 

Jo Jo Jo 

[ de f dx f dyi^l,^ 
Jo Jo 

2tx i-l 







2tt 
1 



dQ / dy 



Jo 



dx 
x=l rr 



dx dx' 



dy 
,*' 

\dy' 



x=0 



dx ' 



and 



d_ 
dy 
1 

27rJo 

r.27r 



5Xo 



9a; 

27r 







-I de I dx I dy- 



"'0 



Jo 

1 



dx 
1 



■0« 



de / dx j dy I dx' j dy' G{x,y,x' ,y')^^ 



Jo 



,,dXo- 
dx 

*:dX^ 
dx' 



1 

2^7o 
1 



de I dx I dyipl^y 



du 
'd^ 



2tt 



-I del dy 



x=0 



dru.,y 

dx 



where we used the following abbreviations: 

X^^Xo{x',y',e), u' = u{x',y',e), u*' ^ u* {x' ,y' ,e), 
and also defined the following functions: 



Jo 



d 



'^*u,x{x,y,e) = I dx I dy' G{x' ,y' ,x,y) — 



C(2;,2/,e) = dx' dy'G{x',y',x,y) — 



u*ix',y',e) 
u*{x',y',e) 



dXo{x',y',e) 
dy' 

dXo{x',y',e) 



dx' 



(A9) 



(AlO) 



(All) 

(A12) 

(A13) 
(A14) 



Here, we note that Eqs. ((29))(|30|) can be derived by applying the Laplacian to Eqs. (jA13p(|A14l) . respectively. In this 
way, the adjoint operator £*{x, y, 0), defined in Eq. (1^ . is obtained as 



£*(x,,,eK(x,,,e) = ^ + ^ + - 



dy 



d_ 

dy 



dx 
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In addition, the adjoint boundary conditions are given by 



du*{x,y, 6) 



dx 



du*{x,y,e) 



x=0 



dx 



= 0, 



x=l 



u*{x,y,@) =u*{x,y,e) 



= 0, 



(A16) 
(A17) 



which represent the Neumann zero boundary condition on x and the Dirichlct zero boundary condition on 
y. In fact, under these adjoint boundary conditions, the bilinear concomitant S[u*(x,y,&),u{x,y,Q)] = 
lu*{x, y, 9), C{x, y, <d)u{x, y, 6)] — y, &)u*{x, y, 0), y, 9)] becomes zero, i.e.. 



u*{x,y,Q),u{x,y, 9) 



27r Jo Jo 
27r Jo Jo 



dy 



dy 



1 
1 

2^ 



de / dx 



de / dx 



dy 



u* — 
dx 

du* ■ 
dx 

* du 

du* 
dy 



x=l 



x=0 
x=l 



x=0 
y=l 



+ 



27r io io 



dx 



y=o 

-I x = l 



2n 



de / dy 



dx 



dXp 
dy 



1- r.9 r 

27r Jo Jo 
27r Jo Jo 



dXo , 



dx 



x=0 

y=o 

j/=0 



x=l 



x=0 



1 

n 

2^ 



d9 / dy 



x=0 
n x=l 



dx 







t dy 


* 

u u 


'0 





:r=0 
e=27r 

= 0. 

e=o 



(A18) 



Each term of the bilinear concomitant S[u*{x, y, 9), u{x, y, 9)] vanishes for the following reasons: the first and second 
terms become zero owing to the Neumann zero boundary condition on x for u and u*, respectively; the third and 

fourth terms, the Dirichlct zero boundary condition on y for u* and u, respectively; the fifth to tenth terms, the 
Dirichlet zero boundary condition on both x and y for tpo, tjJu, i'u,x' ^^'^ '^u,y'^ term, the 27r-periodicity with 

respect to 9 for both u and u*. 
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Hn ® 

FIG. 1; (Color online) (a) Limit-cycle orbit projected on the H11-H22 plane, (b) Wave forms of iJii(O) and -^22(0). The 
Rayleigh number is Ra = 480, and then the collective frequency is ~ 622. 
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FIG. 2: (Color) Snapshots oiTo{x,y,e), i/>o(a;, y, 9), Xo{x,y,e), Uo{x,y,e), and Z{x,y,&) for 6 = 0, tt/2, tt, 3tt/2. 
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FIG. 3: (Color online) Typical shapes of both Xo{x,y,Q) and Z{x,y,Q) with respect to 9 at {x,y) = (0.9,0.9) [TR], (0.1,0.1) 
[BL], (0.1,0.9) [TL], (0.9,0.1) [BR]. 
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FIG. 4: (Color online) (a) Anti-symmetric component of phase coupling function, i.e., ra(0) — r(0) — r(— 0). (b) Time 
evolution of collective phase difference, i.e., \Oi — 62], with the coupling parameter e — 0.05. 



